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We study a voltage biased Josephson junction coupled to two resonators of incommensurate frequencies. 
Using a density approach to analyze the cavity fields and an input-output description to analyze the emitted 
photonic fluxes and their correlation functions, we have shown, both for infinite and finite bandwidth detectors, 
that the emitted radiation is non-classical in the sense that the correlators violates Cauchy-Schwarz inequalities. 

We have also studied the time dependence of the photonic correlations and showed that their line-width becomes 
narrower with the increase of the emission rate approaching from below the threshold limit. 


I. INTRODUCTION 

In atomic cavity quantum electrodynamics (cQED), the 
electric dipole moment of an isolated atom interacts with the 
vacuum state electric field of the cavity. The quantum na¬ 
ture of the electromagnetic field can give rise to coherent 
Rabi oscillations between the atom and the cavity provided 
the relaxation and decoherence rate are smaller than the Rabi 
frequenc)®. Similar physics can be observed when the atom 
is replaced by a two-level solid atate system, often called an 
artificial atom, capacitively coupled to a single mode of the 
microwave cavitjEEl Such paradigmatic systems can be de¬ 
scribed by a Hamiltonian H - //sys + f^cavity + ^int where i/sys 
denotes the Hamiltonian of the atom or artificial atom. The 
interaction of light and matter is encoded in //im while //cavity 
describes the cavity Hamiltonian. 

When the system under study is characterized by a finite 
number of degrees of freedom, the physics is rather well un¬ 
derstood. However, replacing the atom by a mesoscopic con¬ 
ductor with open Fermi reservoirs leads to a dramatically dif¬ 
ferent physics as the whole Fermi sea is affected by the cou¬ 
pling to the cavity. Such setups composed of a quantum con¬ 
ductor coupled to a cavity have been realized experimentally, 
both in experiments using metallic tunnel junction as well 

as in more recent experiments with high-Q microwave cavities 
coupled to either quantum dot or carbon nanotube^^EI] 

For an electric conductor, the dimensionless parameter that 
encodes the light-matter interaction is given by the ratio be¬ 
tween the impedance of the environment and the quantum of 
resistance /?k = h/e^. It can therefore be controlled and sig¬ 
nificantly increased by modulating the impedance of the envi¬ 
ronment. Under such conditions, one can potentially analyze 
the two sides of our system, namely the electronic one by stan¬ 
dard transport measurements or the optical one by measuring 
the power emitted in the cavity. 

The influence of an electromagnetic environment on the 
electronic transport properties of a coherent conductor has 
been studied experimentally rather extensively in the past 
twenty five years. Such a phenomenon is well-understood 
for quantum conductors in the tunneli ng reg ime and is called 
dynamical Coulomb Blockade (DCBj.G^HSlThis occurs when 
a quantum coherent conductor, such as a tunnel junction, is 
inserted into a circuit. When an electron tunnels through a 
junction, this entails sudden voltage variations at the edge 
of the environment which excite its electromagnetic modes. 
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FIG. 1. The sketch of the Josephson junction system coupled to 
two resonators. The Li( 2 ) and C|( 2 ) elements stand for the left 
(right) inductance and capacitance, respectively, and V is the volt¬ 
age biased applied over the JJ. The oscillators frequency is such that 
2eV' = n<i>i + %(ji> 2 . Also, b„ut,i(2) stand for the emitted radiation by the 
left (right) oscillator and collected by the left (right) detector. This ra¬ 
diation is expected to show strong quantum-mechanical correlations 
(see text). 


The backaction of the circuit then affects the charge trans¬ 
fer through the junction and its conductance is reduced in a 
nonlinear way.^S In this respect, the microwave cavity acts 
as structured environment that can absorb photons (or emit 
at finite temperature). DCB has been experimentally probed 
both for non-resonant environment^^® as well as for envi¬ 
ronments formed by resonators,excellent theo¬ 
retical agreement. 

For a Josephson junction (JJ), DCB effects are more dra¬ 
matic since at bias voltages smaller than 2A/e, with A 
the superconducting gap, quasi-particles cannot be excited. 
Therefore the only possible channel to dissipate the energy 
of tunnelling Cooper pairs is to transform it in to photons 
emitted in the electromagnetic environment.®^ In a recent 
experiment,® Hofheinz et al. have observed and character¬ 
ized the radiation associated to the flow of Cooper pair con¬ 
nected to a microwave resonator. In particular, it was evi¬ 
denced two photon processes, for which a single Cooper pair 
tunnelling through the junction emits two photons into the 
environment.® Such two-photons processes which are trig¬ 
gered by tunneling through a quantum conductor are very in¬ 
teresting because they can carry non-classical correlations and 
may therefore offer a source of pairs of correlated photons 
which could be very useful for quantum communications.® 
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The microwave radiation emitted by a Josephson junction 
coupled to a single resonator has been theoretically studied 
recently in various regimes in particular, it has been 
shown by Leppagankas et al^ that the radiation emitted by 
the two photon-processes is non-classical, and that it violates 
a classical Cauchy-Schwarz inequality for two mode power 
cross-correlated fluctuations. Here we theoretically extend 
these considerations by analyzing the radiation emitted by a 
Josephson junction coupled to two resonators of different cav¬ 
ities u)i + 0 ) 2 - The device we study is sketched in Fig. 1. 
We analyze the cotTelations between the photons emitted by 
the JJ in the two resonators at bias matching the condition 
2eV = ha>i + Tia> 2 , but also the power cotTelations of the mi¬ 
crowave fields in the corresponding transmission lines using 
an input-output formalism. 

The paper is organized as follows. In Sec. [I^we introduce 
the setup and the system Hamiltonian of the Josephson junc¬ 
tion coupled to two LC oscillators. In Sec.j^we introduce the 
reduced density matrix describing the coupled system and cal¬ 
culate various observables, such as the average photon number 
and second order coherence functions. In Sec. EYl we utilize 
the input-output theory to characterize the emitted light by the 
two oscillators, for both infinite and finite bandwidths of the 
detectors. We also calculate the time-dependent second order 
coherence functions and calculate the resulting Fano factor of 
the radiation. Finally, in Sec.|V]we end up with some conclu¬ 
sions and outlook. 


II. MODEL HAMILTONIAN 

In the following we will describe the proposed model and 
the associated Hamiltonian. In Fig. [T] we show a schemat¬ 
ics of the voltage biased Josephson junction (JJ) coupled to 
two distinct LC oscillators which, themselves, are coupled to 
a transmission line that can be used to probe each of oscillator 
independently. The model Hamiltonian describing the entire 
setup, including the external lines (which in Fig. [T] quantify 
the emitted radiation) reads; 


third term represents the coupling between the oscillator and 
the environment reads: 

^ ggia^j -t- agj) , (3) 

where g'^ are the coupling strengths. Finally, the last term 
describes the Josephson junction and it’s capacitive coupling 
to the two oscillators: 

Hj = -Ej cos 7? - 2eiv - ^ , (4) 

1 = 1,2 

where Ej, rj, V, N, and Vres are the Josephson energy, the 
phase bias, the voltage bias, the number of Cooper pairs, and 
the voltage induced by the oscillators, respectively. The latter 
is given by = 2,= 1,2 VLs’ with 

= (5) 

Note that we defined the following conjugate variables for the 
cavities and the JJ, respectively; 

[qi,(/)j] ^2ie6ij, (6) 

(7) 

and we have that Ej - Ejo cos (;7r<I)/<I>o), with E the Joseph¬ 
son energy and <t) the flux threaded through the JJ ((Fq = h/2e 
is the flux quantum), namely it can be controlled by control¬ 
ling the flux 11®^. Next we exclude the Copper pairs number 
from the Hamiltonian by performing a time-dependent unitary 
transformation on the full Hamiltonian, nameljEH 

i/tot = UM(t)H,,,ul(t) + i drUNitwiit) (8) 

with 

t/^(f) = e'("-'r+0i+(S2W^ (9) 

where coj = 2eVlfi. That affects both the resonators and JJ 
Hamiltonians so that they become: 


77tot = + hU) + Hj , (1) 

i=l,2 


being the sum of the resonators Hamiltonians, the environ¬ 
ment (the external transmission lines), their coupling to the 
corresponding resonator, and the Josephson junction Hamil¬ 
tonian, respectively. Specifically, 


2C; 


2Li 


( 2 ) 


is the Hamiltonian of the ith (/ = 1,2) resonator, with 0, and 
qi being the phase and charge operators acting on the capaci¬ 
tance C, and inductance L,, respectiveljl^MSOj The second term 
T^env “ the Hamiltonian of the ith environ¬ 

ment, with ; {aq,i) being the photon creation (annihilation) 
operator for the photons in that environment, with cog j and 
q being the energy and their wave-vector, respectively. The 


2Ci \2ej2Li 

Hj = -Ej cos {ojjt + ^ j , 
1 = 1,2 


( 10 ) 

( 11 ) 


where qi - qi + 2eN is still continuos,although N is discrete, 
so that the Hamiltonian of the resonators becomes identical to 
the original one. We can quantize the excitations in the two 
resonators by writing 

= ^[KAa] - ad, qi = i yl¥i(a] + ad ( 12 ) 

• 'y * 1 " 

with Ki — til2mi(jji, — 2mia>fi, mi — (/i/2e )C,-, and a\ (ad 
are the creation (annihilation) operators for the resulting pho¬ 
tons (and which satisfy [a,, aj] = dij). Using this description, 
the oscillators Hamiltonian becomes: 

H'res = OJia]ai , 


(13) 
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with a)i - 1 / yfCiCi. We note that the parameter k, = 
y/LilCi quantifies the range of the coupling regime: for 
Ki 1 the setup is in the weak coupling limit which, com¬ 
bined with small photonic emission rate (oc Ej), is well de¬ 
scribed by the so called P{E) theorySHHI Por /c, > 1 in¬ 
stead, one can reach the strong coupling regime where non- 
linearities and feed-back effects of the electronic transport on 
the cavity (and vice-versa) become manifest. In this paper we 
will only be concerned with the weak coupling limit, since 
most experiments are carried out in this regime. 

Until now, the description was exact, but the coupling be¬ 
tween the JJ and the resonators is strongly non-linear. To 
make further progress, it is instructive to switch to the inter¬ 
action picture with respect to the two cavities, which pertains 
to write: 


Hj{t) = -Ej cos 


cJjt+Yj 

!= 1,2 


(14) 


which is generally time-dependent. In the following we are 
interested in photonic processes that are resonant with uj, 
namely where the energy of the Cooper pairs gained from the 
bias voltage equals the energy of a given number of quanta in 
the two resonators. Here we are interested in the two-photon 
emission processes, one from each cavity, such that the bias 
voltage is tuned to wy = wi H- W 2 - 

Typically, one performs the so called rotating wave approx¬ 
imation (RWA), which amounts to keeping only the terms in 
the above Hamiltonian that do not depend explicitly on time. 
To identify such terms, we expand the argument of the cosine 
function so that we obtain: 


to the single cavity case. We then obtain for the static part of 
the above Hamiltonian the following expression: 




A A mMjii + 1)! ^ 


.1=1,2 Hi 


4 4 . 




Jx (2 ^jKln\) Ji (2 ^|K2n2) 


yjnx 




(18) 


where Ej - £’yexp[-(/ci h-a’ 2)/2], and Jp{z) is the Bessel 
function of the first kind that has the form 


n 


{zl2f’'^P 
nl(n + p)\ 


(19) 


and : ... : means normal ordering of the operators, i.e. all 
annihilation operators on the right and all the creation ones 
on the left. We mention, however, that we still need to add 
the environments and their coupling to the two oscillators to 
complete the approximate setup. While we should in principle 
perform the same RWA on the external modes only (and on 
their coupling), it turns out that for g‘^ « w, we can safely 
utilize the initial coupling Hamiltonian’. In the following we 
analyze the photon emission by the junction by using both the 
density matrix approach and the input-output approach. We 
will calculate quantities such as the average photon number 
(or fluxes) and the second coherence factors that unravel the 
photonic statistics. 


III. DENSITY MATRIX DESCRIPTION 


Hj{t) = 


nZ''""- 

/=1,2 ni,mi 




x(fl^)"'(a,)'"'e''^"'""'^"'' + h.c.] . 


(15) 


Now we perform the RWA by keeping only the static terms, 
namely we impose the condition n\ -mi - ^2 -'«2 = -1- This 
statement is only valid if the two frequencies are incommen¬ 
surate with wy, namely wy/wi + sjr, with s, r e N [and thus 
wy/wy 4^ r/(s - r)]. Otherwise, we would have the condition 


(«! - mi)s ^ (ny - m2){r - ■*) _ 


or 


n\ — mi 


(n2 - m 2 ) - (n2 - m 2 -H 1) 


r 

s] 


(16) 


(17) 


which implies that (^2 - m 2 + l)r/s - q eN. For example, if 
r/s - 1/2, we have that 112 - m 2 + I = 2q and ni-mi - q—\. 
We will not consider this commensurate limit any further here, 
and we leave it as a subject for future analysis. We mention 
that the single-photon and two-photon emission from a single 
cavity have been already discussed theoreticalljl^MSO] 
plemented experimental!}!^. Thus, even if the energy of these 
modes is the same, the resulting processes must not be similar 


In this section we analyze the photons in the cavity by re¬ 
sorting to the density matrix approach, similar to that used in 
Ref Such a description allows us to analyze the photonic 
fields inside each of the two cavities. Let us start by writing 
down the equation of motion for the total photonic density ma¬ 
trix, including the external transmission lines (environments): 

^ ■ ( 20 ) 

1 = 1,2 

so that the reduced density matrix ps(f) describing only the 
oscillators and the JJ can be found by tracing the environment 
degrees of freedom, namely ps{f) — Ji^:env\ptot{t)}. Within the 
Markov approximation, and up to second order in the coupling 
g'^, one obtains the following (Lindblad) form for the system 
density matrix: 

^ = -Cps = + Yj T/ [2a,PsaJ - {«■«/,Psl] , 

1 = 1,2 

which describes the dynamics of the electromagnetic modes 
in the cavity, where y,- = 2jiYjq\g‘q?' is the /th cavity decay 
rate. Note that the expectation values of an operator X reads 
{X{t)) - Trjj,j[Xps(f)], where the trace is taken over the os¬ 
cillators degrees of freedom. In order to characterize the pho¬ 
tonic emission statistics, we define the following second-order 
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correlation functiorl^ 


G®(f,T) = <; np{f)nj{t + t) :) = <: n;(0)e^’'[ps(f)np(0)] :), 

( 21 ) 

which describe the probability to detect a photon in branch p 
at time t and one photon in branch j at time f + t, with r being 
the time delay between the photonic counts. For p j (p - j) 
they correspond the cross-correlation (auto-correlation) of the 
photonic emission, namely detection of photons from different 
(the same) emitter. In the stationary limit, for t oa, the 
correlators depend only on the time delay One defines a 
normalized second correlation function (for t oa) as 


( 2 ) 

a (t t) = - 

<n^(f))<n,(f + T))’ 


( 22 ) 


which, in the stationary limit and for t - 0 (zero time delay) 
becomes: 


^ J ^ V(np,nj)-6pj{nj) 


(23) 


Because the Hamiltonian is not quadratic in the field 

operators, the system of equations for does not close, and we 
find: 


[itj, - Ej -sJiEj : ^ajaj - alaj j 

X [7o(2 ^jKjrij) + 72(2 ^/Kjr^j) 


7i(2y^jnj) 


(28) 


and a similar, but more complicated result for the 
[npnj,Hj^^] commutator (not shown). To make progress, in 
the following we recall again the weak coupling assumption, 
ki ^2 ^ 1, as in the experiments carried out in Ref. [21] and 
consider also that the two cavities are populated on average 
by only a few photons, so that A-,{n,) « 1. The latter con¬ 
dition depend on the parameter range and has to be checked 
self-consistently at the end of the calculation. Under these 
assumptions, we can approximate the system Hamiltonian as 
follows: 


with V(np, rij) - {upfij) - {np){nj) being the variance of the 
field. By accessing the gpj(O) one can infer the degree of 
“quantumness” of the emitted light. Let us now analyze in de¬ 
tail this statement, by comparing the possible outcomes from 
both the classical and quantum descriptions. The correlation 
between the photon emission processes from the two oscilla¬ 
tors is encoded in the variance Var(ni - 112 ) = (((Jni - 5n2)^) > 
0. Assuming {n\) - { 112 } = (n) (equally populated cavities), 
this implies that (dwj) + {dn^} > 2{6n\6n2). For a classical 
field, we can define - {6np6nj}lrP', without normal 

ordering since n, are not operators, but random variables. By 
using the above condition, we obtain for classical states the 
following Cauchy-Schwarz inequalities: 


(2) ^fn(0) + ^®2(0) 


*c,12 


( 0 )< ^ 


^fll(0)^S2(0)- 

(24) 


For a quantum field instead, {dnp5nj) - {: Supdnj :) + 6pj{np) 
and, since in this case ggp/O) - {'■ dnpdnj :), we obtain that 

{6np6nj)l{n)^ = gf^pjiO) + (ll{n})6pj. Thus, for a quantum 
field we obtain the less stringent condition: 


(2) 8%i0) + g%(0) 1 


(25) 


which implies that the Cauchy-Schwarz inequalities can be 
violated. 

Next we analyze explicitly these Cauchy-Schwarz inequal¬ 
ities for our two-photon setup. For that, we first calculate var¬ 
ious average quantities in the stationary limit d{X)ldt - 0, so 
that for the average photon number X - hj and photon number 
product X - hphj, we obtain the following relations: 


fjRWA ^ ■\JkiK2 : {a[a2 + 0102 ) :, ( 29 ) 

neglecting terms of the order OiKini)^. This simple quadratic 
Hamiltonian describes the well known non-degenerate para¬ 
metric amplifieJ^. For the photon number expectation values 
in the two cavities we get: 

{rij) = y {{Hj) + {Up) + 1 ) , (30) 


where p - (2Ej yfKii^ly) for the case of two identical cavities 
(ri = 72 ), so that <ni) = <n 2 ) = <n), with 


(n> = 


2(1-y82)- 


(31) 


We see that /? = 1 is the threshold for an instability at which 
the cavities become coherently populated instead of incoher¬ 
ently, as it happens below this threshold. We will not discuss 
any further the behavior above this threshold, as we are in¬ 
terested here by the incoherent regime. However, note that for 
the quadratic approximation to be valid, the condition *•/(«,) <Si 
must be met, and thus the threshold cannot be achieved within 
this limit as strong quantum fluctuation may change the crit¬ 
ical conditiorl^. Under all these assumptions and after some 
lengthy, but straightforward calculations we obtain: 

^ ^ P^[{n^) + {l+P^){n} + 2{nin2)+pV2] 

- 

+ ^(2<„>+l)+i<„>5^.^, (32) 


i 

'^yj 


inp)^-2^{[np,HT 


(26) 


(27) 


(n^) = A{n^) + B{n\n2) + C , 
{n\n2) - A{n^) + B{n\n2) + D , 


{UpUj) = - 


(33) 

(34) 
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with 


A = 

C = 

D = 




4-^62’ 


>8" 




(1-/?2)(4-/?2) 2(1 -yS2)’ 

y8" , P^ 

(1-;62)(4_^2) 4(1 -yS2)- 


(35) 

(36) 


and Z? = 2A. By further manipulating the above expression, 
we finally obtain: 


{n^)^2{{n)f + {n). 


(nini) = 

which in turn leads to 


<n)(l+4<n» 


(37) 

(38) 


g®(0) = gf^iO) = 2, 


(39) 


( 2 ) 

We mention that for g.. (0) = 2 the radiation resembles the 
thermal radiation, although the system is at T = 0, and it is 
characteristic to the parametric non-degenerate oscillator. The 
cross-correlation instead, gi 2 (®) ~ l/2{n) for (n) 1, mean¬ 

ing strong bunching of the emitted radiation at low emission 
rates oc (n). In order to check the non-classical character of 
the emitted pairs, we calculate the so called noise reduction 
factor {NRF) given by 


NRF = 


Var(ni 


2<n) 


(J2) . 

2 


( 2 ) 

4 ■ 


■2gf^) + l, (40) 


which is NRF > 1 for classical light, when the photon emis¬ 
sion is uncorrelated, and can be 0 < NRF < 1 when quantum 
correlations are manifest. For the few photon limit considered 
here, we find NRF - 1/2, which means the photon pair emis¬ 
sion shows strong quantum correlations, which should be eas¬ 
ily detected in experiments. The very same relation implies 
that the classical Cauchy-Schwarz inequalities in Eq. 24 are 
violated. 

We mention that for larger emission rates (cx Ej) or larger 
a:, ’s, one can go beyond the few photon limit and access higher 
non-linearities of the JJ. Such a regime, although interesting, 
is left for a future stud)® 


IV. INPUT-OUTPUT DESCRIPTION 

While the density matrix approach allows us to calculate 
the field in the cavity, it does not tell us the properties of the 
field exiting in the transmission lines and measured by the de¬ 
tectors. In order to address this issue, we resort to the input- 
output description of the system: an input field is sent to the 
combined system formed by the JJ and the two oscillators, and 
the output field is measured. The input field could even be just 
the quantum vacuum, as it will be assumed in the following. 


In Fig. [T] we show a schematics of the input-output helds. 
We assume that each oscillator is coupled to its external trans¬ 
mission line (environment) by only one side, so that there is 
only one input and output fields, respectively. The relations 
between the input, output, and the cavity fields for oscillator 
a — 1,2 are as followPM 

aa{f} = - ^[OaiO, “ yfy^bin^ait) , (41) 


where ^^(0 [b^uSf)] are the input (output) fields defined as 
follows: 


bin^aif) ~ 


bout,a(f) ~ 


1 

'\j2.7rp(y 

1 

'sj'2.7iP(y 




(43) 

(44) 


with fo < f < fi, and pa - Yjq - <^q,a) is the environment 
a density of states. We mention that for deriving these ex¬ 
pression we assumed the RWA to hold, i.e. we neglected the 
terms of the form a^a^ a (counter rotating terms). 

Because of the second term in Eq. ( |4T] l (the commutator), the 
equation of motion is highly non-linear, and reads: 




X : 


i i ^ 2(2 yjKafla^ ^^7o(2 ^KaTla) 


J\{2 'sjt^afla^ 




(45) 


where in the last line we presuppose the limit Katia 1. 
Moreover, we assume the two oscillators have large Q-factors 
so that their overlap in frequency is negligible, i.e. |tui - W 2 I » 
ri,2- 

Utilizing the same quadratic approximation of the JJ Hamil¬ 
tonian, we arrive at the following equation of motion for the 
cavity held operators: 

aa = -icOaOait) - iEj (0^^''“'^“^^' 

- Yaff(0 - ypY^bin^ait ), (46) 

which implies that the left and right helds are coupled and we 
need to solve the dynamics for both of them at the same time. 
It is instructive to switch to the Fourier space in order to obtain 
the needed relationships 

1 r“ 

aa{t) - —— I d(jje“^‘aa(cLi ), (47) 

J-oc 

and the commutation relations: 


WoioS), aa(cj')] - 6{aj - a >'), (48) 


and similarly for the other operators. This in turn allows us 
to write the following equation relating the input and cavity 
helds: 

l(0J -OJa)- Y 

= ^l%bin,a{(^) . (49) 


aa(u) -I- Ej y[i^a\(u\ +U2- co) 
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Manipulating this relation and its hermitian conjugate, and 
utilizing the relation between the input, output and cavity 
fields in Eq. (|44l), we obtain: 


a„(w) = Aa{cS)bin,a{o^) + (wi + W2 - w) , (50) 

+ yl^Ba{aj)bl^ -{a)i + W2 - w), (51) 

with 

Aaiuj) =- ^ -, (52) 

-1{U> -OJa) + Y+ 

iEj yjKaKaJa 

•/ \ To- E^KaKo 

-U)o,) + ^ + i(„_^^)_^./2 

(53) 

We are now in position to use the above findings to calculate 
the relevant observables and correlators. However, we will 
analyze both the case where the detectors have infinite band¬ 
width, namely all photons emitted are collected, as well as 
the finite bandwidth case, and state the differences compared 
to the density matrix approach. Let us start with the infinite 
bandwidth case, and then discuss briefly the implications of 
finite bandwidth detection. 


[((w - w„) - 


A. Infinite bandwidth detection 

Here we assume that the efficiency of the detectors is unity 
and that their bandwidth is infinite, thus they are collecting all 
the emitted photons. The outgoing photonic flux from oscilla¬ 
tor a reads: 


From these expression, we can readily evaluate the zero time 
delay (t = 0) second-order coherence functions: 


(bLAt)bLjt)b 

out,a (0^OHf,O'(0) 

ibout,a(Obout,a(t)){bg^fi^(t)bout^a.{t)} 


2 + 


2 , 

yfnj 

V2r’ 


a' - a 
a' - a. 


(57) 


We see that in the infinite bandwidth limit we obtain again 
the same result we found from the density matrix approach 
for the autocorrelation function, g® = 2, while gj^^O) cx l/T 
for r/y « 1. Note that experimentally it is measured not the 
cavity field, but the external photon flux and the above results 
are due to the well-known relationships 


Ta = {blut,a(t)bout,a(t)} = 7a{al(t)aa(t )), (58) 

{blu,J0blu,a(t)bout,a(t)bout,a(0) = rl(al(t)al(t)aa(t)aa(t)}. 

(59) 


which hold in the infinite bandwidth case. Once again, devi¬ 
ations from the quadratic Hamiltonian lead to changes in the 
second order coherence factors, of the order (kY/ y). Since we 
consider the few photon regime, such corrections are negligi¬ 
ble. It worth mentioning that the due to the correspondence 
in Eq. ( |59] l the Cauchy-Schwarz inequalities are in the same 
way as described in the previous section, thus the emitted light 
from the cavities is non-classical. 

Next we calculate the time-dependence of the g„„,(T) func¬ 
tion in order to reveal the dynamics of the two resonators, and 
which is defined as follows: 

^(2) ^ 

■ (60) 


r«(0 - {blu,Jt)bo,,,a{t)) = — du\ 

^71 ^_oo %J—oo 

'X<.bl^ta(^)bout,oW)) = r daj\Ba(aj)\^ = Ir, 

Y2nJ-oc 2(1-yS^) 

(54) 


which, by assuming yi - y 2 = y implies Eq, s E, while the 
autocorrelation and cross-correlation functions at t = 0 re¬ 
spectively read: 


G?i(0) = <^L,.(0^L,„(f)^„„r..(f)^o„,,<.(f)> 

,,2 f /^oo i2 „q4 


% 

n 


f 


du\Ba{u))\ 


nir 


(l-^2)2’ 


(55) 


gS(0) = 

= — ^ dco J" dco'[B^(a))A'^(a)i + a>2 - a>) 


X Ba(oj')Acy(a>i + a»2 - oj') + |-6q'(w)| \Ba{co')\ 


'\\ 2 -\ 


ny^pW +0^) 
2{l-p2p 


( 56 ) 


The auto-correlation and cross-correlation second order co¬ 
herence functions are given respectively by: 


Gal(r) = <^Lr,a(f + T)bl^,^g(f)bout,a{t)boutAt + T)> 

Try^jB^ r. -.-/ . .. . _ . , .v 2 


2(1 - p 2)2 


1 H- e (cosh(£ jt) + p sinh(£yx)) 


- {b',u,.ait + r)bl,,,(t)bou,At)bou,At + T)} 
ny^P^ 


2(1 -p^P 


I + e [p cosh(£ jt) + sinh(£yT)^ 


( 61 ) 


which leads to the following normalized coherence functions: 

[^cosh(jSy|T|)-H sinh(ySy|T|)]^ , 


g'^),(T) = 1 -H 1 

o aa’ ' ' 


e '>''4 [cosh(jSy|T|) -i-/3sinh(ySy|T|)]2 


, a - a. 
(62) 


These functions are witnesses of the dynamics of of oscil- 
laltors as they reveal their effective linewidth and thus their 
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FIG. 2. Left plot: The time dependence of the second-order auto¬ 
correlation function gaa(T) for different values of J3 oc Ej. Right 
plot: The frequency dependence of the auto-correlation noise func¬ 
tion Saa(T) for different values of p oc Ej. In both plots, the red, blue, 
brown, green, and yellow curves corresponds tOyS = 0.1,0.3,0.7,0.8, 
and 0.95, respectively. Here the time scales are all expressed in terms 
of 7 , and we recall that T oc 


FIG. 3. Left plot: The time dependence of the second-order cross¬ 
correlation function g®(T) for different values of p oc Ej. Right 
plot: The frequency dependence of the auto-correlation noise func¬ 
tion Saa{T) for different values of p oc Ej. In both plots, the red, blue, 
brown, green, and yellow curves corresponds to/3 = 0.1,0.3,0.7,0.8, 
and 0.95, respectively. 


trend towards the instability point >6=1. Having found these 
functions, we can directly relate them to the non-symmetrized 
photonic frequency noise, defined as follows: 


S a a '— 


£/T<ra(T)ra,.(0))e 


r 

X oo 

oo 


(63) 


with Fq-Ct) = bl^i,^^(T)bout,a(T) representing the photonic flux 
operator out of the oscillator a. That allows us to extract the 
noise to signal ratio, or the Fano factor Fao' = Saa'{(o>)/T. 
This offers an estimate for the number of photons that are cor¬ 
related in the emission proces^, and it is given by: 

X OO 

• (64) 

OO 

Using the expression for the photonic rate T in Eq. (HI, and 
expressing the coefficient p in terms of this quantity, we obtain 
for the auto-correlation Fano factor: 


F 


aa' 


5F 4F2 

7 r , 
5F 4F2 

2 + - + -T- , 


r 


r 


a' - a 
a' - a'. 


(65) 


We see that for F/y » 1 (but so that /cF/y < 1), F oc F^, 
which is a signature for strong photon bunching. We note that 
for a Poisson process F - I, while for a thermal distribution 
F = 2. 

In Fig. 1^ we plot the dependence of the functions g® (r) 
on the time delay (left) and the autocorrelated noise Saa((^) 
(right). We see that for r —> oo 1, while from the 

right plot we see that the linewidth becomes narrower as the 
emission rate oc Ej increases to values close to the instability 
threshold. In Fig. [3] on the other hand, we plot the depen- 
dence of the functions on the time delay (left) and the 

cross-correlated noise SaaC^o) (right), which shows a similar 
behavior. 


B. Finite bandwidth detection 


In this section we discuss the properties of the outgoing 
light when the detection is taken within a finite frequency win¬ 
dow Aw around both wi and W 2 . For simplicity, we assume 
the same frequency window for both oscillators. We already 
expect that if Aw » y the results to be practically the same as 
in the previous section. However, when Aw < y, one expects 
the photon statistics to be affected, and we aim at quantifying 
such changes. 

The photon flux detected within the bandwidth Aw reads: 


^^out,a 

out,aify) /\(jj 

jP 

arctan 


I-P 


y 

yfh: 

arctan 


f 

Oixi, 




ijjQ~l!:^col2 

(it-)! 


l+y6 


t/w|F„(w)|^ 

( 66 ) 


where r - KojIj. The rate becomes Eq. ( |5^ for Aw » y up 
to corrections {jjKijSf’. Eor r < 1, we obtain: 


^ O' ~ 


4 j 62 Aw 

(l->62)2 


+ 0[r^]. 


(67) 


Next we calculate the second-order correlation functions. 
We obtain: 


G^AO) = — 


arctan 


(l^) arctan 


1-/6 


l+>6 


G%i0) = 


y2>62 


arctan () arctan () 
( 1+^)2 + ( 1_^)2 


( 68 ) 

(69) 


(2^ 

which results again in g„„(0) = 2, independent on the band¬ 
width Aw, while for the cross-correlation coefficient we ob¬ 
tain: 


slriO) = 1 + 


[(1 -Hj6)arctan(Y^)H-(l -/?) arctan 
[(1 -H >6) arctan (y^)-( 1 -/?) arctan 

( 70 ) 






























8 




FIG. 4. Left: The cross-correlation factor ^(0) as a function of the 
bandwidth r = Aw/y for various values of the rate F. The red, blue, 
brown, green, yellow curves correspond to F = 0.1, 0.3, 0.5, 1, and 
3, respectively. Right: The bi 2 function in Eq. l |71[ (, together with the 
approximate expression in Eq. Here all curves, for which we 
used the same parameters as in the left plot, lie on the same line, and 
are fitted well by the approximate result. 



For completeness, let us also investigate the time- 
dependence of the second-order correlation functions for a fi¬ 
nite bandwidth detection. They are defined as in Eq. ( |60l l, with 
the second-order correlation functions and rate being calcu¬ 
lated over a finite bandwidth. While even in this case there are 
possible analytical expressions for the ^(t) functions, they 
are too lengthy and uninspiring to be displayed. However, in 
the left (right) plot in Fig. we show galrCr) as a 

function of the delay time t for various values of r. We see 
that the finite bandwidth strongly modifies the decay of the 
correlation functions for r < 1, which now shows correlations 
over a longer time scale, of the order 1 / Aw. Note that the ef¬ 
fective linewidth of the oscillators in the presence of the JJ is 
of the order of y(l - j3^), which needs to be compared with 
1/Aw. The cross-correlation Fano factor is also modified, and 
for 1 it reads T'j 2 = 2 + 5r/(r), with 




n[33r + 40r^ -H 15r^ -H 15(1 -H arctan r] 
30(1 + r^)[r H- (1 + r'^) arctan r]^ 


(73) 


This function behaves as f(r) ^ InjSr {f(r) « 1 H- 8/37rr^) 
for r —> 0 (r ^ oa). 


FIG. 5. Left (right): The auto-correlation (cross-correlation) factor 
g£m,r(r) ,.(t)] as a function of the time delay t for different values 
of r = Aw/y. The black, red, blue, brown, green, yellow curves 
correspond to r = oo, 3, 1, 0.5, and 0.2, respectively, and we choose 
F = 0.1. 


which instead dependents on the bandwidth Aw. Let us ana¬ 
lyze this expression in more detail. In experiments, one mea¬ 
sures the photonic rate exiting the device, instead of the bare 
emission rate p. Thus, we should calculate gW(0) as a func¬ 
tion of Fr^a- Assuming the emission rate is such that the sys¬ 
tem is far below the instability threshold, p <Si 1, we found the 
general form: 


( 2 ) _ , ^'2 
^Ur(O) -ai2+ • 

1 ry 


(71) 


with fli 2 a! 2 for jS < 1 and all values of r, while bn is given 
by 


?12 



arctan r - 


r 

2 + 


(72) 


In Fig. I^we show the dependence of g® ^(0) (left) and bn 
(right) on r for different values of the measured rate F. We 
see that bn increases monotonically with increasing Aw, and 
that it saturates to b‘^™ - ylnjl. To conclude the zero-time 
delay discussion, we see that the input-output results match 
the our findings from the density matrix approach, but that 
considering a finite bandwidth affects the correlation function 
which shows less bunching (i.e. less correlated emission) as 
the bandwidth is decreased to smaller values than the natural 
bandwidths of the oscillators F. 


V. DISCUSSIONS AND CONCLUSIONS 

In conclusion, we have studied the radiation emitted by a 
voltage biased JJ into two LC oscillators when the the bias 
is set such that 2eV - tia>i + tia> 2 . We have employed both 
the density matrix approach to study the cavity fields, and 
the input-output description to analyze the emitted photonic 
fluxes in the so called weak coupling regime characterized by 
an environmental impedance much smaller than Rk- Specif¬ 
ically, we have calculated both the photon number and the 
photonic correlations (second-order coherence function) and 
showed, by proving that a Cauchy-Schwarz inequality is vi¬ 
olated, that the emitted radiation is non-classical. In order 
to analyze the effective dynamics of the oscillators, we also 
calculated the time-dependence of the photonic correlations 
and showed that the their linewidth becomes narrower as the 
emission rate increases, signaling the approach to the thresh¬ 
old instability. We also briefly considered the effect of a finite 
bandwidth detection and showed that the qualitative features 
stay the same, but that the violation of the Cauchy-Schwarz in¬ 
equality becomes less pronounced. In the future, it would be 
interesting to address also the strong coupling regime together 
with the region around the threshold instability, as well as the 
full counting statistics of the emitted photons, as discussed in 
Ref Ell 
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Note added. During the final completion of this paper, we arXiv: 1503.00545 using the density matrix approach for the 

became aware of similar results obtained by A. Armour et al. cavity fields. 
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